Data were loaded from spreadsheets obtained by Utah DWQ:
# tributary data
trib <- read.csv("./Data/Processed/ul_Tribdata_wqp_processed_2020-12-17.csv")
# facility data
fac <- read.csv("./Data/Processed/ul_Facilitydata_wqp_processed.csv")
# watershed spatial data
wshds <- read_sf("./Data/GIS/ulSubWshdsWGS84.shp")
# sampling activity logs
activitylog <- read.csv("./Data/Processed/UL_Trib_SampleActivities_exported 12-17-20.csv")
# groundwater data
gw <- read.csv("./Data/Processed/ul_wqp_gwdata_processed_2021-01-07.csv")
# precipitation and evaporation data
p_et <- read.xlsx("./Data/Processed/ULEFDCFlowBoundaries_WY2006-2018.xlsx",
sheet = "Precip-ET", startRow = 2)
# DMR data
dmr <- read.xlsx("./Data/Processed/DMR_Analysis_v3_updated 12-7-20.xlsx",
sheet = "DMR.Formatted")
dmr_LUT <- read.xlsx("./Data/Processed/DMR_Analysis_v3_updated 12-7-20.xlsx",
sheet = "LUT_SiteList", startRow = 2, cols = 1:5)
Watersheds were divided into monitored and unmonitored and arranged in clockwise order.
wshds_monitored <- wshds %>%
filter(Monitored_ == "YES")
wshds_monitored$WS_Name <-
factor(wshds_monitored$WS_Name,
levels = c("Lake Mountain - Israel Canyon", "Dry Creek - Saratoga",
"Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough", "Currant Creek"))
levels(wshds_monitored$WS_Name)
## [1] "Lake Mountain - Israel Canyon" "Dry Creek - Saratoga"
## [3] "Lehi Spring Creek" "American Fork River"
## [5] "Timp SSD" "Lindon Drain"
## [7] "Powell Slough Major" "Provo River"
## [9] "Mill Race" "Hobble Creek"
## [11] "Dry Creek - Spanish Fork" "Spanish Fork River"
## [13] "4000 South Drain Spanish Fork" "Benjamin Slough"
## [15] "Currant Creek"
Tributary dataset was filtered for:
Flow data were converted to cubic feet/sec (cfs) when units originally appeared as gallon/min, cubic meters/sec, and million gallons/day.
trib$ActivityStartDate <- as.Date(trib$ActivityStartDate)
# unique(trib$CharacteristicName)
trib$ResultMeasureValue <- as.numeric(as.character(trib$ResultMeasureValue))
## Warning: NAs introduced by coercion
trib_wq <- trib %>%
# Set nondetects (qualifier code U) as 1/2 MDL.
# Qualifier code J is above detection but below reporting. Values retained as-is.
# REVISIT: What does qualifier code Q mean? One sample for ammonium at 0.062
mutate(ResultMeasureValue = case_when(MeasureQualifierCode == "U" ~ MethodDetectionLevel/2,
TRUE ~ ResultMeasureValue)) %>%
filter(CharacteristicName %in%
c("Ammonia-nitrogen as N", "Carbon", "Nitrate as N", "Organic carbon",
"Ammonium", "Carbonate", "Inorganic carbon", "Nitrate", "Nitrogen",
"Phosphorus", "Ammonia and ammonium",
"Inorganic nitrogen (nitrate and nitrite)", "Nitrite", "Orthophosphate",
"Kjeldahl nitrogen", "Nitrogen, mixed forms (NH3), (NH4), organic, (NO2) and (NO3)",
"Organic Nitrogen", "Ammonia-nitrogen", "Phosphate-phosphorus",
"Total Kjeldahl nitrogen")) %>%
filter(ResultMeasure.MeasureUnitCode %in% c("ueq/L", "mg/l", "mg/l as N", "<NA>", "mg/l as P")) %>%
select(MonitoringLocationIdentifier:ResultSampleFractionText,
OrganizationIdentifier:MonitoringLocationTypeName, LatitudeMeasure,
LongitudeMeasure, ResultMeasureValue, ResultMeasure.MeasureUnitCode,
ResultStatusIdentifier) %>%
mutate(Month = month(ActivityStartDate),
Year = year(ActivityStartDate)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
filter(!is.na(ResultMeasureValue)) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultMeasure.MeasureUnitCode) %>%
unite("Variable", CharacteristicName:ResultSampleFractionText)
trib_flow <- trib %>%
filter(CharacteristicName %in% c("Flow", "Stream flow, instantaneous", "Stream flow, mean. daily")) %>%
droplevels() %>%
#filter(ResultStatusIdentifier != "Provisional") %>% # can filter out provisional data
select(MonitoringLocationIdentifier:ResultSampleFractionText,
OrganizationIdentifier:MonitoringLocationTypeName, LatitudeMeasure,
LongitudeMeasure, ResultMeasureValue, ResultMeasure.MeasureUnitCode,
ResultStatusIdentifier) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultSampleFractionText,
-MonitoringLocationTypeName, -CharacteristicName) %>%
group_by(MonitoringLocationIdentifier, ActivityStartDate, OrganizationIdentifier,
OrganizationFormalName, MonitoringLocationName, LatitudeMeasure,
LongitudeMeasure, ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
summarise(ResultMeasureValue = mean(ResultMeasureValue)) %>%
mutate(Month = month(ActivityStartDate),
Year = year(ActivityStartDate)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
pivot_wider(names_from = ResultMeasure.MeasureUnitCode, values_from = ResultMeasureValue) %>%
# USGS data has both cfs and cms, always duplicated. remove cfs.
select(-'ft3/s')
## `summarise()` has grouped output by 'MonitoringLocationIdentifier', 'ActivityStartDate', 'OrganizationIdentifier', 'OrganizationFormalName', 'MonitoringLocationName', 'LatitudeMeasure', 'LongitudeMeasure', 'ResultMeasure.MeasureUnitCode'. You can override using the `.groups` argument.
colnames(trib_flow)[11:14] <- c("m3.sec", "cfs", "mgd", "g.min")
# convert everything to cfs
trib_flow$mgd <- trib_flow$mgd/0.646317
trib_flow$g.min <- trib_flow$g.min*0.002228017
trib_flow$m3.sec <- trib_flow$m3.sec*35.3147
trib_flow_long <- trib_flow %>%
pivot_longer(cols = c("cfs", "mgd", "g.min", "m3.sec"),
names_to = "orig.units", values_to = "Flow.cfs") %>%
drop_na(Flow.cfs)
# conversions seem to match up
ggplot(trib_flow_long, aes(x = ActivityStartDate, y = Flow.cfs, color = orig.units)) +
geom_point(alpha = 0.5) +
scale_y_log10() +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.8)
## Warning: Transformation introduced infinite values in continuous y-axis
A sample type of “11” indicates a site was not sampled. A sample type of “10” indicates there was no flow. The latter should appear in the flow data as a zero; this was verified by inspecting the two datasets.
Locations and dates when sites were not sampled:
## [1] "LOWER S FK PROVO R AT GAGING STATIION"
## [2] "N FK PROVO R AB CNFL / PROVO R AT WILDWOOD"
## [3] "CURRANT CK AB MONA RES"
## [4] "Timpanogos Effluent below constructed duck ponds"
## [5] "Dry Ck Near Utah Lake-WLA"
## [6] "DRY CK EAST TRIBUTARY"
## [7] "AMERICAN FK CK 2.5MI S OF AM FK CITY"
## [8] "Powell Slough WMA North Outfall to Utah Lake"
## [9] "Powell Slough WMA South Outfall to Utah Lake"
## [10] "Saratoga Springs at Cedar Valley"
## [1] "2011-01-26" "2014-02-19" "2012-01-19" "2012-02-15" "2012-03-21"
## [6] "2012-04-10" "2012-06-20" "2012-07-25" "2011-10-27" "2011-11-16"
## [11] "2011-12-06" "2011-12-07" "2017-05-17" "2017-07-10" "2017-12-20"
## [16] "2017-12-22" "2018-01-10" "2018-02-15" "2018-03-13" "2018-04-19"
## [21] "2018-05-17" "2018-11-20" "2019-05-14" "2019-06-13" "2019-06-17"
## [26] "2019-07-08" "2019-08-12" "2019-09-24" "2019-09-23" "2019-10-08"
## [31] "2019-11-04" "2019-12-10"
Tributary dataset was filtered for:
Flow data were converted to cubic feet/sec (cfs) when units originally appeared as gallon/min, cubic meters/sec, and million gallons/day.
fac$datetime <- as.POSIXct(fac$datetime)
#unique(fac$CharacteristicName)
fac_wq <- fac %>%
# Set nondetects (qualifier code U) as 1/2 MDL.
# Qualifier code J is above detection but below reporting. Values retained as-is.
# REVISIT: What does qualifier code Q mean? One sample for ammonium at 0.062
mutate(ResultMeasureValue = case_when(MeasureQualifierCode == "U" ~ MethodDetectionLevel/2,
TRUE ~ ResultMeasureValue)) %>%
filter(CharacteristicName %in%
c("Ammonia-nitrogen as N", "Carbon", "Nitrate as N", "Organic carbon",
"Ammonium", "Carbonate", "Inorganic carbon", "Nitrate", "Nitrogen",
"Phosphorus", "Ammonia and ammonium",
"Inorganic nitrogen (nitrate and nitrite)", "Nitrite", "Orthophosphate",
"Kjeldahl nitrogen", "Nitrogen, mixed forms (NH3), (NH4), organic, (NO2) and (NO3)",
"Organic Nitrogen", "Ammonia-nitrogen", "Phosphate-phosphorus",
"Total Kjeldahl nitrogen")) %>%
filter(ResultMeasure.MeasureUnitCode %in% c("ueq/L", "mg/l", "mg/l as N", "<NA>", "mg/l as P")) %>%
select(MonitoringLocationIdentifier:LongitudeMeasure, ResultMeasureValue,
ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
mutate(Month = month(datetime),
Year = year(datetime)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
filter(!is.na(ResultMeasureValue)) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultMeasure.MeasureUnitCode) %>%
unite("Variable", CharacteristicName:ResultSampleFractionText)
# match date field with trib
colnames(fac_wq)[2] <- "ActivityStartDate"
fac_flow <- fac %>%
filter(CharacteristicName == "Flow") %>%
#c("Height, gage", "Stream flow, instantaneous", "Stream flow, mean. daily")) %>%
droplevels() %>%
#filter(ResultStatusIdentifier != "Provisional") %>% # can filter out provisional data
select(MonitoringLocationIdentifier:LongitudeMeasure, ResultMeasureValue,
ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
select(-SampleDepthValue, -RelativeDepth, -ResultSampleFractionText,
-MonitoringLocationTypeName, -CharacteristicName) %>%
group_by(MonitoringLocationIdentifier, datetime, OrganizationIdentifier,
OrganizationFormalName, MonitoringLocationName, LatitudeMeasure,
LongitudeMeasure, ResultMeasure.MeasureUnitCode, ResultStatusIdentifier) %>%
summarise(ResultMeasureValue = mean(ResultMeasureValue)) %>%
mutate(Month = month(datetime),
Year = year(datetime)) %>%
filter(Year >= 2010) %>%
droplevels() %>%
pivot_wider(names_from = ResultMeasure.MeasureUnitCode, values_from = ResultMeasureValue)
## `summarise()` has grouped output by 'MonitoringLocationIdentifier', 'datetime', 'OrganizationIdentifier', 'OrganizationFormalName', 'MonitoringLocationName', 'LatitudeMeasure', 'LongitudeMeasure', 'ResultMeasure.MeasureUnitCode'. You can override using the `.groups` argument.
fac_flow$datetime <- as.Date(fac_flow$datetime)
colnames(fac_flow)[2] <- "ActivityStartDate"
colnames(fac_flow)[11:14] <- c("mgd", "cfs", "g.min", "m3.sec")
# convert everything to cfs
fac_flow$mgd <- fac_flow$mgd/0.646317
fac_flow$g.min <- fac_flow$g.min*0.002228017
fac_flow$m3.sec <- fac_flow$m3.sec*35.3147
fac_flow_long <- fac_flow %>%
pivot_longer(cols = c("cfs", "mgd", "g.min", "m3.sec"),
names_to = "orig.units", values_to = "Flow.cfs") %>%
drop_na(Flow.cfs)
Tributary and facility water quality and flow datasets were combined.
USGS_flow <- readNWISdv(siteNumbers = c("10146400", "10150500",
"10153100", "10163000"),
parameterCd = "00060",
startDate = "2010-01-01")
USGS_flow_processed <- USGS_flow %>%
mutate(WS_Name = case_when(site_no == "10146400" ~ "Currant Creek",
site_no == "10150500" ~ "Spanish Fork River",
site_no == "10153100" ~ "Hobble Creek",
site_no == "10163000" ~ "Provo River")) %>%
select(WS_Name, Date, X_00060_00003)
colnames(USGS_flow_processed) <- c("WS_Name", "ActivityStartDate", "Flow.USGS.cfs")
Sites located downstream of diversion
Sites that may be backwatered at times
Intermittent sites
Effluent-dominated sites
Sites included in ULWQS SAP
trib_fac_sites <- trib_fac_wq_flow %>%
ungroup() %>%
select(MonitoringLocationIdentifier, MonitoringLocationName,
MonitoringLocationTypeName, LatitudeMeasure, LongitudeMeasure) %>%
distinct() %>%
mutate(SiteType = ifelse(MonitoringLocationTypeName %in%
c("Canal Drainage", "Canal Irrigation", "Canal Transport"), "Canal",
ifelse(MonitoringLocationTypeName %in% c("River/Stream", "Stream", "Stream: Ditch"), "Stream",
ifelse(MonitoringLocationTypeName %in% c("Facility Industrial", "Facility Other"), "Facility",
"Storm Sewer"))))
trib_fac_sites_sf <- st_as_sf(trib_fac_sites, coords = c("LongitudeMeasure", "LatitudeMeasure"),
crs = 4326, dim = "XY") #CRS = WGS84
trib_fac_sites_intersect <- trib_fac_sites_sf %>%
st_intersection(wshds) %>%
left_join(trib_fac_wq_flow, .,)
unique(trib_fac_sites_intersect$MonitoringLocationName[is.na(trib_fac_sites_intersect$WS_Name)])
## [1] "Provo River Delta Restoration"
## [2] "JORDAN R AT UTAH L OUTLET U121 XING"
## [3] "Powell Slough WMA North Outfall to Utah Lake"
## [4] "Powell Slough WMA South Outfall to Utah Lake"
## [5] "Dry Ck Near Utah Lake-WLA"
## [6] "Unnamed flow at 4000 West @ mouth"
## [7] "American Fork River at mouth"
## [8] "Timpanogos WWTP at Utah Lake mouth"
## [9] "Lindon Drain at Utah Lake Inlet"
## [10] "Beer Creek at Utah Lake mouth"
## [11] "Hobble Creek at Utah Lake mouth"
## [12] "Mill Race at mouth 1.25 miles west of I-15"
## [13] "Dry Creek North Fork @ Confluence"
## [14] "Dry Creek South Fork @ Confluence"
unique(trib_fac_sites_intersect$MonitoringLocationIdentifier[is.na(trib_fac_sites_intersect$WS_Name)])
## [1] "UTAHDWQ_WQX-4917343" "UTAHDWQ_WQX-4994790" "UTAHDWQ_WQX-4995210"
## [4] "UTAHDWQ_WQX-4995230" "UTAHDWQ_WQX-4996040" "WFWQC_UT-4917712"
## [7] "WFWQC_UT-4994958" "WFWQC_UT-4995043" "WFWQC_UT-4995075"
## [10] "WFWQC_UT-4995210" "WFWQC_UT-4995467" "WFWQC_UT-4996040"
## [13] "WFWQC_UT-4996096" "WFWQC_UT-4996536" "WFWQC_UT-4996047"
## [16] "WFWQC_UT-4996048"
unique(trib_fac_sites_intersect$WS_Name)
## [1] NA "Lehi Drain 2"
## [3] "Lake Mountain - Israel Canyon" "Dry Creek - Saratoga"
## [5] "Lehi Spring Creek" "American Fork River"
## [7] "Timp SSD" "Lindon Drain"
## [9] "Vineyard Drain 6" "Powell Slough Major"
## [11] "Currant Creek" "Benjamin Slough"
## [13] "Spanish Fork River" "Dry Creek - Spanish Fork"
## [15] "Provo Bay 3" "Provo Bay 5"
## [17] "Hobble Creek" "Mill Race"
## [19] "Provo River" "4000 South Drain Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Powell Slough WMA North Outfall to Utah Lake"] <- "Powell Slough Major"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Powell Slough WMA South Outfall to Utah Lake"] <- "Powell Slough Major"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Dry Ck Near Utah Lake-WLA"] <- "Dry Creek - Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Mill Race at mouth 1.25 miles west of I-15"] <- "Mill Race"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Hobble Creek at Utah Lake mouth"] <- "Hobble Creek"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Beer Creek at Utah Lake mouth"] <- "Benjamin Slough"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Lindon Drain at Utah Lake Inlet"] <- "Lindon Drain"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Timpanogos WWTP at Utah Lake mouth"] <- "Timp SSD"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"American Fork River at mouth"] <- "American Fork River"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Unnamed flow at 4000 West @ mouth"] <- "4000 South Drain Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"JORDAN R AT UTAH L OUTLET U121 XING"] <- "Jordan River"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"DRY CK EAST TRIBUTARY"] <- "Dry Creek - Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"Drainage Canal 0.5 mile bl I-15 at about 2500 West, Springville"] <- "Dry Creek - Spanish Fork"
trib_fac_sites_intersect$WS_Name[trib_fac_sites_intersect$MonitoringLocationName ==
"PROVO STATION 6-WLA"] <- "Mill Race"
trib_fac_sites_wsds <- trib_fac_sites_intersect %>%
ungroup() %>%
select(MonitoringLocationIdentifier, OrganizationIdentifier:LongitudeMeasure,
SiteType:Shape_Area) %>%
distinct() %>%
arrange(WS_Name)
# write.csv(trib_fac_sites_wsds, file = "./Data/Processed/trib_fac_sites_initial.csv",
# row.names = FALSE)
trib_fac_sites_counts <- trib_fac_sites_intersect %>%
ungroup() %>%
group_by(MonitoringLocationIdentifier, OrganizationIdentifier, OrganizationFormalName,
MonitoringLocationName, MonitoringLocationTypeName, LatitudeMeasure,
LongitudeMeasure, SiteType, OBJECTID_1, OBJECTID, TARGET_FID,
Type, wsa_sqkm, WS_Name, Monitored_, WS_NAME_1, WS_NAME_12,
Shape_Leng, Shape_Le_1, Shape_Area) %>%
tally() %>%
ungroup() %>%
select(MonitoringLocationIdentifier, n)
# trib_fac_sites_SD <- read.csv("./Data/Processed/trib_fac_sites.csv")
#
# trib_fac_sites_combined <- full_join(trib_fac_sites_SD, trib_fac_sites_counts)
# write.csv(trib_fac_sites_combined, file = "./Data/Processed/trib_fac_sites_counts.csv",
# row.names = FALSE)
trib_fac_sites_wsds <- trib_fac_sites_wsds %>%
full_join(., trib_fac_sites_counts)
# write.csv(trib_fac_sites_wsds, file = "./Data/Processed/trib_fac_sites_counts.csv",
# row.names = FALSE)
trib_fac_sites_intersect_downstream <- trib_fac_sites_intersect %>%
filter(MonitoringLocationIdentifier %in%
c("UTAHDWQ_WQX-4995312", "UTAHDWQ_WQX-4995310", # Currant Creek, duplicate sites at same lat/lon
"UTAHDWQ_WQX-4995465", "WFWQC_UT-4995467", # Benjamin Slough
"UTAHDWQ_WQX-5919910", "WFWQC_UT-4917712", # 4000 South Drain Spanish Fork
"UTAHDWQ_WQX-4995578", "WFWQC_UT-4995575", # Spanish Fork River
"UTAHDWQ_WQX-4996040", "WFWQC_UT-4996040", # Dry Creek - Spanish Fork
"UTAHDWQ_WQX-4996042", "UTAHDWQ_WQX-4996044", # Dry Creek - Spanish Fork
"WFWQC_UT-4996100", "UTAHDWQ_WQX-4996100", "WFWQC_UT-4996096",
"UTAHDWQ_WQX-4996275", "WFWQC_UT-4996275", # Hobble Creek
"WFWQC_UT-4996536", "WFWQC_UT-4996540", "UTAHDWQ_WQX-4996540", "UTAHDWQ_WQX-4996566", # Mill Race
"WFWQC_UT-4996680", "UTAHDWQ_WQX-4996680", # Provo River
"UTAHDWQ_WQX-4995230", "WFWQC_UT-4995210", "UTAHDWQ_WQX-4995210", # Powell Slough Major
"WFWQC_UT-4995075", "UTAHDWQ_WQX-4995120", # Lindon Drain
"WFWQC_UT-4995043", "UTAHDWQ_WQX-4995041", "UTAHDWQ_WQX-4995038", # Timp SSD
"WFWQC_UT-4994958", "UTAHDWQ_WQX-4994960", # American Fork River
"WFWQC_UT-4994948", "UTAHDWQ_WQX-4994950", #Lehi Spring Creek
"UTAHDWQ_WQX-4994804", # Dry Creek - Saratoga
"UTAHDWQ_WQX-4994792")) %>% #Lake Mountain - Israel Canyon
mutate(Backwater = ifelse(MonitoringLocationIdentifier %in%
c("WFWQC_UT-4917712", "WFWQC_UT-4994958", "UTAHDWQ_WQX-4995465",
"WFWQC_UT-4995467", "WFWQC_UT-4996536", "UTAHDWQ_WQX-4995210",
"UTAHDWQ_WQX-4995230", "WFWQC_UT-4995210", "WFWQC_UT-4995575"),
"Yes", "No")) %>%
left_join(., USGS_flow_processed)
unique(trib_fac_sites_intersect_downstream$WS_Name)
## [1] "Lake Mountain - Israel Canyon" "Dry Creek - Saratoga"
## [3] "Lehi Spring Creek" "American Fork River"
## [5] "Timp SSD" "Lindon Drain"
## [7] "Powell Slough Major" "Currant Creek"
## [9] "Benjamin Slough" "Spanish Fork River"
## [11] "Dry Creek - Spanish Fork" "Hobble Creek"
## [13] "Mill Race" "Provo River"
## [15] "4000 South Drain Spanish Fork"
# check that flow at USGS gage is in line with measured flow at chem sites, when data exist together.
# Some cases in Spanish Fork River where measured flow was near zero but USGS flow was not. Otherwise lines up on 1:1 line.
# Will proceed to replace missing flow with USGS flow in those
# ggplot(subset(trib_fac_sites_intersect_downstream, WS_Name %in%
# c("Hobble Creek", "Spanish Fork River", "Provo River", "Currant Creek")),
# aes(x = Flow.cfs, y = Flow.USGS.cfs, color = WS_Name)) +
# geom_point() +
# geom_abline(slope = 1, intercept = 0)
trib_fac_sites_intersect_downstream$Flow.cfs[is.na(trib_fac_sites_intersect_downstream$Flow.cfs)] <-
trib_fac_sites_intersect_downstream$Flow.USGS.cfs[is.na(trib_fac_sites_intersect_downstream$Flow.cfs)]
# remove zero flow samples from watersheds that do not go dry
# then, calculate daily loads
trib_fac_sites_intersect_downstream <- trib_fac_sites_intersect_downstream %>%
mutate(Flow.cfs = case_when(WS_Name == "Lehi Spring Creek" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Timp SSD" & OrganizationIdentifier == "WFWQC_UT" &
Flow.cfs == 0 ~ NA_real_,
WS_Name == "Powell Slough Major" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Provo River" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Mill Race" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Hobble Creek" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Dry Creek - Spanish Fork" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Spanish Fork River" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "4000 South Drain Spanish Fork" & Flow.cfs == 0 ~ NA_real_,
WS_Name == "Benjamin Slough" & Flow.cfs == 0 ~ NA_real_,
TRUE ~ Flow.cfs)) %>%
mutate(load.TN.kgd = TN * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TDN.kgd = TDN * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TP.kgd = TP * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TDP.kgd = TDP * Flow.cfs * 28.3168 * 86400 / 1000000,
load.TOC.kgd = TOC * Flow.cfs * 28.3168 * 86400 / 1000000,
load.DOC.kgd = DOC * Flow.cfs * 28.3168 * 86400 / 1000000)
trib_fac_sites_intersect_downstream$load.TN.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TDN.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TP.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TDP.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.TOC.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$load.DOC.kgd[trib_fac_sites_intersect_downstream$Flow.cfs == 0] <- 0
trib_fac_sites_intersect_downstream$SiteType[is.na(trib_fac_sites_intersect_downstream$SiteType)] <- "Not specified"
trib_fac_sites_downstream <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
select(MonitoringLocationIdentifier, MonitoringLocationName, OrganizationIdentifier,
MonitoringLocationTypeName, LatitudeMeasure, LongitudeMeasure, WS_Name) %>%
distinct()
trib_fac_sites_downstream_sf <- st_as_sf(trib_fac_sites_downstream, coords = c("LongitudeMeasure", "LatitudeMeasure"),
crs = 4326, dim = "XY") #CRS = WGS84
trib_fac_sites_downstream_sf$WS_Name <-
factor(trib_fac_sites_downstream_sf$WS_Name,
levels = c("Lake Mountain - Israel Canyon", "Dry Creek - Saratoga",
"Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough", "Currant Creek"))
trib_fac_sites_intersect_downstream$WS_Name <-
factor(trib_fac_sites_intersect_downstream$WS_Name,
levels = c("Lake Mountain - Israel Canyon", "Dry Creek - Saratoga",
"Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough", "Currant Creek"))
# Create two datasets: one starting in 2010, one starting in 2015. The later period of record is the main dataset.
trib_fac_sites_intersect_downstream_early <- trib_fac_sites_intersect_downstream
trib_fac_sites_intersect_downstream <- trib_fac_sites_intersect_downstream %>%
filter(Year >= 2015)
REVISIT
TO DO: Use area vs. elevation relationship to establish average monthly lake area, with variability quantified.
REVISIT: P and ET only available 2005-2018 from EFDC output. We can only use four years of data if we restrict to 2015 onward. Should we use all years?
dmr_load <- dmr %>%
filter(perm_feature_type_desc == "External Outfall") %>%
filter(Nuts.ParmName %in%
c("Nitrogen, ammonia total [as N]", "Nitrogen, Kjeldahl, total [as N]",
"Nitrogen, nitrate total [as N]", "Phosphate, ortho [as P]",
"Phosphorus, total [as P]")) %>%
select(npdes_id, Year, Month,
value_type_desc, dmr_value_nmbr, dmr_unit_desc, statistical_base_short_desc,
NPES_Name, NPDES.Type, NPDESID.Date, Nuts.ParmName, Load) %>%
na.omit(Load) %>%
group_by(npdes_id, Year, Month,
value_type_desc, dmr_value_nmbr, dmr_unit_desc, statistical_base_short_desc,
NPES_Name, NPDES.Type, NPDESID.Date, Nuts.ParmName) %>%
summarise(Load = mean(Load)) %>%
ungroup() %>%
separate(col = NPDESID.Date, into = c(NA, "Date")) %>%
pivot_wider(names_from = "Nuts.ParmName", values_from = "Load")
## `summarise()` has grouped output by 'npdes_id', 'Year', 'Month', 'value_type_desc', 'dmr_value_nmbr', 'dmr_unit_desc', 'statistical_base_short_desc', 'NPES_Name', 'NPDES.Type', 'NPDESID.Date'. You can override using the `.groups` argument.
dmr_load$Date <- as.Date(as.numeric(dmr_load$Date), origin = "1899-12-30")
colnames(dmr_load) <- c("NPDES.ID", "Year", "Month", "value.type", "value",
"units", "stat", "NPDES.Name", "Type", "Date",
"TP", "NH3", "PO4", "TKN", "NO3")
# # 30 day avg most available for TP
# ggplot(dmr_load, aes(x = NPDES.Name, y = TP, fill = stat)) +
# geom_boxplot()
#
# # 30 day avg most available for PO4
# ggplot(dmr_load, aes(x = NPDES.Name, y = PO4, fill = stat)) +
# geom_boxplot()
#
# # 30 day avg most available for TKN
# ggplot(dmr_load, aes(x = NPDES.Name, y = TKN, fill = stat)) +
# geom_boxplot()
#
# # data very limited for NO3. may want to use TKN instead.
# ggplot(dmr_load, aes(x = NPDES.Name, y = NO3, fill = stat)) +
# geom_boxplot()
#
# # several stats available for NH3. may want to use TKN instead.
# ggplot(dmr_load, aes(x = NPDES.Name, y = NH3, fill = stat)) +
# geom_boxplot()
dmr_load <- dmr_load %>%
filter(stat == "30DA AVG") %>%
# load is in ton/yr. convert to metric tons/mo
mutate(TP.tonmo = TP/(12*1.10231),
PO4.tonmo = PO4/(12*1.10231),
TKN.tonmo = TKN/(12*1.10231),
NH3.tonmo = NH3/(12*1.10231),
NO3.tonmo = NO3/(12*1.10231)) %>%
select(-TP, -PO4, -TKN, -NH3, -NO3) %>%
mutate(MonitoringLocationIdentifier =
case_when(NPDES.Name == "MCWANE DUCTILE-UTAH" ~ "UTAHDWQ_WQX-4996430",
NPDES.Name == "OREM CITY CORP" ~ "UTAHDWQ_WQX-4995250",
NPDES.Name == "PacifiCorp Energy" ~ "UTAHDWQ_WQX-4995124",
NPDES.Name == "PAYSON POWER PROJECT" ~ "UTAHDWQ_WQX-4995480",
NPDES.Name == "PROVO CITY CORP" ~ "UTAHDWQ_WQX-4996560",
NPDES.Name == "SALEM CITY CORP" ~ "UTAHDWQ_WQX-4995440",
NPDES.Name == "SPANISH FORK CITY CORP" ~ "UTAHDWQ_WQX-4996020",
NPDES.Name == "SPRINGVILLE- CITY OF" ~ "UTAHDWQ_WQX-4996280",
NPDES.Name == "TIMPANOGOS SPECIAL SERVICE DIS" ~ "UTAHDWQ_WQX-4995040"))
dmr_sites <- dmr_load %>%
select(NPDES.Name, MonitoringLocationIdentifier) %>%
distinct() %>%
left_join(., trib_fac_sites) %>%
filter(!is.na(MonitoringLocationTypeName))
## Joining, by = "MonitoringLocationIdentifier"
dmr_sites_sf <- st_as_sf(dmr_sites, coords = c("LongitudeMeasure", "LatitudeMeasure"),
crs = 4326, dim = "XY") #CRS = WGS84
dmr_sites_sf <- st_intersection(dmr_sites_sf, wshds)
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## Warning: attribute variables are assumed to be spatially constant throughout all
## geometries
The top graph depicts all tributary and facility monitoring sites. The second graph depicts each monitored watershed individually from the NW side of the lake moving clockwise, with downstream monitoring sites marked. The third graph depicts all watersheds with only downstream monitoring sites marked.
ggplot(dmr_load, aes(x = Date, y = TP.tonmo)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(NPDES.Name), nrow = 2) +
scale_y_log10() +
labs(y = "TP load (metric ton/mo)")
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Removed 897 rows containing missing values (geom_point).
ggplot(dmr_load, aes(x = Date, y = PO4.tonmo)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(NPDES.Name), nrow = 2) +
scale_y_log10() +
labs(y = "TDP load (metric ton/mo)")
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Removed 1247 rows containing missing values (geom_point).
ggplot(dmr_load, aes(x = Date, y = TKN.tonmo)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(NPDES.Name), nrow = 2) +
scale_y_log10() +
labs(y = "TKN load (metric ton/mo)")
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Removed 1015 rows containing missing values (geom_point).
ggplot(dmr_load, aes(x = as.factor(Month), y = TP.tonmo)) +
geom_boxplot() +
facet_wrap(vars(NPDES.Name), nrow = 2) +
scale_y_log10() +
labs(y = "TP load (metric ton/mo)", x = "Month")
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Removed 909 rows containing non-finite values (stat_boxplot).
ggplot(dmr_load, aes(x = as.factor(Month), y = PO4.tonmo)) +
geom_boxplot() +
facet_wrap(vars(NPDES.Name), nrow = 2) +
scale_y_log10() +
labs(y = "TDP load (metric ton/mo)", x = "Month")
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Removed 1259 rows containing non-finite values (stat_boxplot).
ggplot(dmr_load, aes(x = as.factor(Month), y = TKN.tonmo)) +
geom_boxplot() +
facet_wrap(vars(NPDES.Name), nrow = 2) +
scale_y_log10() +
labs(y = "TKN load (metric ton/mo)", x = "Month")
## Warning: Transformation introduced infinite values in continuous y-axis
## Warning: Removed 1029 rows containing non-finite values (stat_boxplot).
There is only one inflow monitored in Lake Mountain - Israel Canyon. Estimates are based on this site.
Flow
Nutrient loads
There is only one inflow monitored in Dry Creek - Saratoga. Estimates are based on this site.
Flow
Nutrient loads
There are two sites in Lehi Spring Creek. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
There are two sites in American Fork River. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Flow
Nutrient loading
REVISIT East side tributary (UTAHDWQ_WQX-4995041) has measurable concentrations but no flow (one near-zero flow, one flow duplicated from 4995038). Scott’s recommendation was to average concentrations at 4995041 and 4995038 and apply them to the flow at 4995038, but this assumes that the volume-weighted average is the same as the arithmetic average (i.e., that the volumes at both sites are the same). I’m not sure if this assumption is valid. For now, I will ignore the east side tributary and calculate loads based on the sites that have flow data.
There are two sites in Lindon Drain. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Flow
Nutrient loading
Powell Slough has both North and South outfalls. Need to add together load from North and South to compute total load. WFWQC samples only at the north outfall.
There are two site codes in Provo River (one DWQ and one WFWQC), but these are in the same location. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Flow
Nutrient loading
4996536 represents the sum total of loading from the watershed. 4996540 and 4996566 should be added together to represent total flow.
Flow
Nutrient loading
4996096 and 4996100 are equivalent sites. 4996275 is a separate flow representing Springville and represents its own load.
Flow
Nutrient loading
4996040 is the primary site in this watershed and represents the sum total of loads. 4996042 and 4996044 represent equivalent sites for the East tributary but do not include the south fork. East tributary sites will be removed and the primary site used for loads.
There are two sites in Spanish Fork River. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
There are two sites in 4000 South Drain Spanish Fork. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
There are two sites in Benjamin Slough. DWQ site is upstream of WFWQC site. Will compare loads to determine if the estimates can be combined.
Flow
Nutrient loading
Two sites, UTAHDWQ_WQX-4995310 and UTAHDWQ_WQX-4995312, are located at the same coordinates. These sites are pooled in load calculations, as they represent the same location.
Flow
Nutrient loads
Monitoring locations
# Lake Mountain - Israel Canyon
loading_LakeMountainIsraelCanyon <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Lake Mountain - Israel Canyon") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Dry Creek - Saratoga
loading_DryCreekSaratoga <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Dry Creek - Saratoga") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Lehi Spring Creek
loading_LehiSpringCreek <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Lehi Spring Creek") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# American Fork River
loading_AmericanForkRiver <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "American Fork River") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Timp SSD
loading_TimpSSD <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Timp SSD") %>%
filter(MonitoringLocationIdentifier != "UTAHDWQ_WQX-4995041") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Lindon Drain
loading_LindonDrain <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Lindon Drain") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Powell Slough
loading_PowellSloughMajor_initial <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Powell Slough Major") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# calculate load from N and S outfalls with WFWQC data from N outfall (DWQ site only at S outfall)
loading_PowellSloughMajor_WFWQCloads <- loading_PowellSloughMajor_initial %>%
filter(MonitoringLocationIdentifier %in% c("UTAHDWQ_WQX-4995230", "WFWQC_UT-4995210")) %>%
group_by(WS_Name, Month) %>%
summarise(Flow.cfs = sum(Flow.cfs, na.rm = TRUE),
load.TN.kgd = sum(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = sum(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = sum(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = sum(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = sum(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = sum(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = sum(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonmo = sum(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonmo = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = sum(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonmo = sum(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonmo = sum(load.DOC.tonmo, na.rm = TRUE)) %>%
mutate(MonitoringLocationIdentifier = "WFWQC_UT-4995210 + UTAHDWQ_WQX-4995230",
OrganizationIdentifier = "WFWQC_UT") %>%
select(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month,
Flow.cfs, load.TN.kgd, load.TDN.kgd, load.TP.kgd, load.TDP.kgd,
load.TOC.kgd, load.DOC.kgd, load.TN.tonmo, load.TDN.tonmo, load.TP.tonmo,
load.TDP.tonmo, load.TOC.tonmo, load.DOC.tonmo)
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
# calculate load from N and S outfalls with DWQ data
loading_PowellSloughMajor_UTAHDWQloads <- loading_PowellSloughMajor_initial %>%
filter(MonitoringLocationIdentifier %in% c("UTAHDWQ_WQX-4995230", "UTAHDWQ_WQX-4995210")) %>%
group_by(WS_Name, Month) %>%
summarise(Flow.cfs = sum(Flow.cfs, na.rm = TRUE),
load.TN.kgd = sum(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = sum(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = sum(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = sum(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = sum(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = sum(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = sum(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonmo = sum(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonmo = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = sum(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonmo = sum(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonmo = sum(load.DOC.tonmo, na.rm = TRUE)) %>%
mutate(MonitoringLocationIdentifier = "UTAHDWQ_WQX-4995210 + UTAHDWQ_WQX-4995230",
OrganizationIdentifier = "UTAHDWQ_WQX") %>%
select(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month,
Flow.cfs, load.TN.kgd, load.TDN.kgd, load.TP.kgd, load.TDP.kgd,
load.TOC.kgd, load.DOC.kgd, load.TN.tonmo, load.TDN.tonmo, load.TP.tonmo,
load.TDP.tonmo, load.TOC.tonmo, load.DOC.tonmo)
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
#combine loading from different organizations
loading_PowellSloughMajor <- rbind(loading_PowellSloughMajor_UTAHDWQloads, loading_PowellSloughMajor_WFWQCloads)
# Provo River
loading_ProvoRiver <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Provo River") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Mill Race
loading_MillRace_initial <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Mill Race") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# calculate mean of 6536 and 6540
loading_MillRace_WFWQCloads <- loading_MillRace_initial %>%
filter(MonitoringLocationIdentifier %in% c("WFWQC_UT-4996536", "WFWQC_UT-4996540")) %>%
group_by(WS_Name, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = mean(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonmo = mean(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonmo = mean(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonmo = mean(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonmo = mean(load.DOC.tonmo, na.rm = TRUE)) %>%
mutate(MonitoringLocationIdentifier = "WFWQC_UT-4996536 + WFWQC_UT-4996540",
OrganizationIdentifier = "WFWQC_UT") %>%
select(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month,
Flow.cfs, load.TN.kgd, load.TDN.kgd, load.TP.kgd, load.TDP.kgd,
load.TOC.kgd, load.DOC.kgd, load.TN.tonmo, load.TDN.tonmo, load.TP.tonmo,
load.TDP.tonmo, load.TOC.tonmo, load.DOC.tonmo)
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
# calculate sum of 6540 and 6566
loading_MillRace_UTAHDWQloads <- loading_MillRace_initial %>%
filter(MonitoringLocationIdentifier %in% c("UTAHDWQ_WQX-4996540", "UTAHDWQ_WQX-4996566")) %>%
group_by(WS_Name, Month) %>%
summarise(Flow.cfs = sum(Flow.cfs, na.rm = TRUE),
load.TN.kgd = sum(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = sum(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = sum(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = sum(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = sum(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = sum(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = sum(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonmo = sum(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonmo = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = sum(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonmo = sum(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonmo = sum(load.DOC.tonmo, na.rm = TRUE)) %>%
mutate(MonitoringLocationIdentifier = "UTAHDWQ_WQX-4996540 + UTAHDWQ_WQX-4996566",
OrganizationIdentifier = "UTAHDWQ_WQX") %>%
select(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month,
Flow.cfs, load.TN.kgd, load.TDN.kgd, load.TP.kgd, load.TDP.kgd,
load.TOC.kgd, load.DOC.kgd, load.TN.tonmo, load.TDN.tonmo, load.TP.tonmo,
load.TDP.tonmo, load.TOC.tonmo, load.DOC.tonmo)
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
#combine loading from different organizations
loading_MillRace<- rbind(loading_MillRace_UTAHDWQloads, loading_MillRace_WFWQCloads)
# Hobble Creek
loading_HobbleCreek_initial <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Hobble Creek") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# calculate load from 6100/6096 and 6275 with WFWQC data
# first, take average of 6100 and 6096 (equivalent sites)
# second, sum 6100/6096 and 6275
loading_HobbleCreek_WFWQCloads_initial <- loading_HobbleCreek_initial %>%
filter(MonitoringLocationIdentifier %in% c("WFWQC_UT-4996096", "WFWQC_UT-4996100")) %>%
group_by(WS_Name, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = mean(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonmo = mean(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonmo = mean(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonmo = mean(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonmo = mean(load.DOC.tonmo, na.rm = TRUE)) %>%
mutate(MonitoringLocationIdentifier = "WFWQC_UT-4996096 + WFWQC_UT-4996100",
OrganizationIdentifier = "WFWQC_UT") %>%
select(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month,
Flow.cfs, load.TN.kgd, load.TDN.kgd, load.TP.kgd, load.TDP.kgd,
load.TOC.kgd, load.DOC.kgd, load.TN.tonmo, load.TDN.tonmo, load.TP.tonmo,
load.TDP.tonmo, load.TOC.tonmo, load.DOC.tonmo)
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
loading_HobbleCreek_WFWQCloads <- loading_HobbleCreek_initial %>%
rbind(., loading_HobbleCreek_WFWQCloads_initial) %>%
filter(MonitoringLocationIdentifier %in% c("WFWQC_UT-4996096 + WFWQC_UT-4996100", "WFWQC_UT-4996275")) %>%
group_by(WS_Name, Month) %>%
summarise(Flow.cfs = sum(Flow.cfs, na.rm = TRUE),
load.TN.kgd = sum(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = sum(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = sum(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = sum(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = sum(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = sum(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = sum(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonmo = sum(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonmo = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = sum(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonmo = sum(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonmo = sum(load.DOC.tonmo, na.rm = TRUE)) %>%
mutate(MonitoringLocationIdentifier = "WFWQC_UT-4996096/4996100 + WFWQC_UT-4996275",
OrganizationIdentifier = "WFWQC_UT") %>%
select(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month,
Flow.cfs, load.TN.kgd, load.TDN.kgd, load.TP.kgd, load.TDP.kgd,
load.TOC.kgd, load.DOC.kgd, load.TN.tonmo, load.TDN.tonmo, load.TP.tonmo,
load.TDP.tonmo, load.TOC.tonmo, load.DOC.tonmo)
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
# calculate load from 6100 and 6275 with DWQ data
loading_HobbleCreek_UTAHDWQloads <- loading_HobbleCreek_initial %>%
filter(MonitoringLocationIdentifier %in% c("UTAHDWQ_WQX-4996100", "UTAHDWQ_WQX-4996275")) %>%
group_by(WS_Name, Month) %>%
summarise(Flow.cfs = sum(Flow.cfs, na.rm = TRUE),
load.TN.kgd = sum(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = sum(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = sum(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = sum(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = sum(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = sum(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = sum(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonmo = sum(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonmo = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = sum(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonmo = sum(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonmo = sum(load.DOC.tonmo, na.rm = TRUE)) %>%
mutate(MonitoringLocationIdentifier = "UTAHDWQ_WQX-4996100 + UTAHDWQ_WQX-4996275",
OrganizationIdentifier = "UTAHDWQ_WQX") %>%
select(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month,
Flow.cfs, load.TN.kgd, load.TDN.kgd, load.TP.kgd, load.TDP.kgd,
load.TOC.kgd, load.DOC.kgd, load.TN.tonmo, load.TDN.tonmo, load.TP.tonmo,
load.TDP.tonmo, load.TOC.tonmo, load.DOC.tonmo)
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
#combine loading from different organizations
loading_HobbleCreek <- rbind(loading_HobbleCreek_UTAHDWQloads, loading_HobbleCreek_WFWQCloads)
# Dry Creek - Spanish Fork
loading_DryCreekSpanishFork <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Dry Creek - Spanish Fork") %>%
filter(MonitoringLocationIdentifier != "UTAHDWQ_WQX-4996042" &
MonitoringLocationIdentifier != "UTAHDWQ_WQX-4996044") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Spanish Fork
loading_SpanishFork <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Spanish Fork River") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# 4000 South Drain Spanish Fork
loading_4000SouthDrain <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "4000 South Drain Spanish Fork") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Benjamin Slough
loading_BenjaminSlough <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Benjamin Slough") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
# Currant Creek
# set site codes to the same site
trib_fac_sites_intersect_downstream$MonitoringLocationIdentifier[trib_fac_sites_intersect_downstream$MonitoringLocationIdentifier == "UTAHDWQ_WQX-4995312"] <-
"UTAHDWQ_WQX-4995310"
loading_CurrantCreek <- trib_fac_sites_intersect_downstream %>%
ungroup() %>%
filter(WS_Name == "Currant Creek") %>%
group_by(WS_Name, OrganizationIdentifier, MonitoringLocationIdentifier, Month) %>%
summarise(Flow.cfs = mean(Flow.cfs, na.rm = TRUE),
load.TN.kgd = mean(load.TN.kgd, na.rm = TRUE),
load.TDN.kgd = mean(load.TDN.kgd, na.rm = TRUE),
load.TP.kgd = mean(load.TP.kgd, na.rm = TRUE),
load.TDP.kgd = mean(load.TDP.kgd, na.rm = TRUE),
load.TOC.kgd = mean(load.TOC.kgd, na.rm = TRUE),
load.DOC.kgd = mean(load.DOC.kgd, na.rm = TRUE),
load.TN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TN.kgd * 30/1000,
Month %in% c(2) ~ load.TN.kgd * 28/1000),
load.TDN.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDN.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDN.kgd * 30/1000,
Month %in% c(2) ~ load.TDN.kgd * 28/1000),
load.TP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TP.kgd * 30/1000,
Month %in% c(2) ~ load.TP.kgd * 28/1000),
load.TDP.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TDP.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TDP.kgd * 30/1000,
Month %in% c(2) ~ load.TDP.kgd * 28/1000),
load.TOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.TOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.TOC.kgd * 30/1000,
Month %in% c(2) ~ load.TOC.kgd * 28/1000),
load.DOC.tonmo = case_when(Month %in% c(1, 3, 5, 7, 8, 10, 12) ~ load.DOC.kgd * 31/1000,
Month %in% c(4, 6, 9, 11) ~ load.DOC.kgd * 30/1000,
Month %in% c(2) ~ load.DOC.kgd * 28/1000)) %>%
distinct()
## `summarise()` has grouped output by 'WS_Name', 'OrganizationIdentifier', 'MonitoringLocationIdentifier', 'Month'. You can override using the `.groups` argument.
DMR Data
# Timp SSD
loading_TimpSSD_dmr <- dmr_load %>%
ungroup() %>%
filter(NPDES.Name == "TIMPANOGOS SPECIAL SERVICE DIS") %>%
group_by(NPDES.Name, Month) %>%
summarise(load.TP.tonmo = mean(TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(PO4.tonmo, na.rm = TRUE),
load.TKN.tonmo = mean(TKN.tonmo, na.rm = TRUE)) %>%
distinct() %>%
mutate(WS_Name = "Timp SSD")
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
# Lindon Drain
# Calculate monthly loading from energy facility
loading_LindonDrain_dmr <- dmr_load %>%
ungroup() %>%
filter(NPDES.Name == "PacifiCorp Energy") %>%
group_by(NPDES.Name, Month) %>%
summarise(load.TP.tonmo = mean(TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(PO4.tonmo, na.rm = TRUE),
load.TKN.tonmo = mean(TKN.tonmo, na.rm = TRUE)) %>%
distinct() %>%
mutate(WS_Name = "Lindon Drain")
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
# Powell Slough
# Calculate monthly loading from Orem WWTP
loading_PowellSloughMajor_dmr <- dmr_load %>%
ungroup() %>%
filter(NPDES.Name == "OREM CITY CORP") %>%
group_by(NPDES.Name, Month) %>%
summarise(load.TP.tonmo = mean(TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(PO4.tonmo, na.rm = TRUE),
load.TKN.tonmo = mean(TKN.tonmo, na.rm = TRUE)) %>%
distinct() %>%
mutate(WS_Name = "Powell Slough Major")
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
# Mill Race
# Calculate monthly loading from Provo WWTP
loading_MillRace_dmr <- dmr_load %>%
ungroup() %>%
filter(NPDES.Name == "PROVO CITY CORP") %>%
group_by(NPDES.Name, Month) %>%
summarise(load.TP.tonmo = mean(TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(PO4.tonmo, na.rm = TRUE),
load.TKN.tonmo = mean(TKN.tonmo, na.rm = TRUE)) %>%
distinct() %>%
mutate(WS_Name = "Mill Race")
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
# Hobble Creek
# Calculate loading from Springville WWTP
loading_HobbleCreek_dmr <- dmr_load %>%
ungroup() %>%
filter(NPDES.Name == "SPRINGVILLE- CITY OF") %>%
group_by(NPDES.Name, Month) %>%
summarise(load.TP.tonmo = mean(TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(PO4.tonmo, na.rm = TRUE),
load.TKN.tonmo = mean(TKN.tonmo, na.rm = TRUE)) %>%
distinct() %>%
mutate(WS_Name = "Hobble Creek")
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
# Dry Creek Spanish Fork
# Calculate loading from Spanish Fork WWTP
loading_DryCreekSpanishFork_dmr <- dmr_load %>%
ungroup() %>%
filter(NPDES.Name == "SPANISH FORK CITY CORP") %>%
group_by(NPDES.Name, Month) %>%
summarise(load.TP.tonmo = mean(TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(PO4.tonmo, na.rm = TRUE),
load.TKN.tonmo = mean(TKN.tonmo, na.rm = TRUE)) %>%
distinct() %>%
mutate(WS_Name = "Dry Creek - Spanish Fork")
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
# Benjamin Slough
# Calculate loading from Payson and Salem. Add loads together.
loading_BenjaminSlough_dmr <- dmr_load %>%
ungroup() %>%
filter(NPDES.Name %in% c("SALEM CITY CORP", "PAYSON POWER PROJECT")) %>%
group_by(NPDES.Name, Month) %>%
summarise(load.TP.tonmo = mean(TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = mean(PO4.tonmo, na.rm = TRUE),
load.TKN.tonmo = mean(TKN.tonmo, na.rm = TRUE)) %>%
distinct() %>%
group_by(Month) %>%
summarise(load.TP.tonmo = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonmo = sum(load.TDP.tonmo, na.rm = TRUE),
load.TKN.tonmo = sum(load.TKN.tonmo, na.rm = TRUE)) %>%
mutate(WS_Name = "Benjamin Slough",
NPDES.Name = "SALEM CITY CORP & PAYSON POWER PROJECT") %>%
select(NPDES.Name, Month, load.TP.tonmo, load.TDP.tonmo, load.TKN.tonmo, WS_Name)
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
loading_AllWatersheds_dmr <-
rbind(loading_BenjaminSlough_dmr, loading_DryCreekSpanishFork_dmr,
loading_HobbleCreek_dmr, loading_LindonDrain_dmr, loading_MillRace_dmr,
loading_PowellSloughMajor_dmr, loading_TimpSSD_dmr)
loading_AllWatersheds_dmr$WS_Name <-
factor(loading_AllWatersheds_dmr$WS_Name,
levels = c("Timp SSD", "Lindon Drain", "Powell Slough Major",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Benjamin Slough"))
loading_AllWatersheds_dmr_annual <- loading_AllWatersheds_dmr %>%
group_by(NPDES.Name, WS_Name) %>%
summarise(load.TP.tonyr = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonyr = sum(load.TDP.tonmo, na.rm = TRUE),
load.TKN.tonyr = sum(load.TKN.tonmo, na.rm = TRUE))
## `summarise()` has grouped output by 'NPDES.Name'. You can override using the `.groups` argument.
loading_AllWatersheds_dmr2_annual <- loading_AllWatersheds_dmr_annual %>%
ungroup() %>%
mutate(OrganizationIdentifier = "DMR",
Flow.cfs = NA,
load.TN.tonyr = load.TKN.tonyr,
load.TDN.tonyr = NA,
load.TOC.tonyr = NA,
load.DOC.tonyr = NA) %>%
select(-NPDES.Name, -load.TKN.tonyr) %>%
select(WS_Name, OrganizationIdentifier, Flow.cfs, load.TN.tonyr, load.TDN.tonyr,
load.TP.tonyr, load.TDP.tonyr, load.TOC.tonyr, load.DOC.tonyr)
# write.csv(loading_AllWatersheds_dmr_annual,
# file = "./Data/Processed/loading_AllWatersheds_dmr_annual.csv",
# row.names = FALSE)
loading_AllWatersheds <- rbind(
loading_LakeMountainIsraelCanyon, loading_DryCreekSaratoga,
loading_LehiSpringCreek, loading_AmericanForkRiver, loading_TimpSSD,
loading_LindonDrain, loading_PowellSloughMajor, loading_ProvoRiver,
loading_MillRace, loading_HobbleCreek, loading_DryCreekSpanishFork,
loading_SpanishFork, loading_4000SouthDrain, loading_BenjaminSlough,
loading_CurrantCreek)
# calculate difference between UTAHDWQ and WFWQC load estimates
loading_AllWatersheds_differences <- loading_AllWatersheds %>%
filter(WS_Name %in% c("Lehi Spring Creek", "American Fork River", "Timp SSD",
"Lindon Drain", "Powell Slough Major", "Provo River",
"Mill Race", "Hobble Creek", "Dry Creek - Spanish Fork",
"Spanish Fork River", "4000 South Drain Spanish Fork",
"Benjamin Slough")) %>%
select(-(load.TN.kgd:load.DOC.kgd)) %>%
group_by(WS_Name, Month) %>%
summarise(load.TN.difference = load.TN.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TN.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TDN.difference = load.TDN.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TDN.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TP.difference = load.TP.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TP.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TDP.difference = load.TDP.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TN.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.TOC.difference = load.TOC.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.TOC.tonmo[OrganizationIdentifier == "WFWQC_UT"],
load.DOC.difference = load.DOC.tonmo[OrganizationIdentifier == "UTAHDWQ_WQX"] -
load.DOC.tonmo[OrganizationIdentifier == "WFWQC_UT"])
## `summarise()` has grouped output by 'WS_Name', 'Month'. You can override using the `.groups` argument.
loading_AllWatersheds_annual <- loading_AllWatersheds %>%
group_by(WS_Name, OrganizationIdentifier) %>%
summarise(Flow.cfs = sum(Flow.cfs),
load.TN.tonyr = sum(load.TN.tonmo, na.rm = TRUE),
load.TDN.tonyr = sum(load.TDN.tonmo, na.rm = TRUE),
load.TP.tonyr = sum(load.TP.tonmo, na.rm = TRUE),
load.TDP.tonyr = sum(load.TDP.tonmo, na.rm = TRUE),
load.TOC.tonyr = sum(load.TOC.tonmo, na.rm = TRUE),
load.DOC.tonyr = sum(load.DOC.tonmo, na.rm = TRUE))
## `summarise()` has grouped output by 'WS_Name'. You can override using the `.groups` argument.
loading_Total_annual <- loading_AllWatersheds_annual %>%
filter(OrganizationIdentifier == "UTAHDWQ_WQX") %>%
ungroup() %>%
summarise(Flow.cfs = sum(Flow.cfs),
load.TN.tonyr = sum(load.TN.tonyr, na.rm = TRUE),
load.TDN.tonyr = sum(load.TDN.tonyr, na.rm = TRUE),
load.TP.tonyr = sum(load.TP.tonyr, na.rm = TRUE),
load.TDP.tonyr = sum(load.TDP.tonyr, na.rm = TRUE),
load.TOC.tonyr = sum(load.TOC.tonyr, na.rm = TRUE),
load.DOC.tonyr = sum(load.DOC.tonyr, na.rm = TRUE))
loading_AllWatersheds_annual_combined <-
rbind(loading_AllWatersheds_annual, loading_AllWatersheds_dmr2_annual)
# write.csv(loading_AllWatersheds_annual, file = "./Data/Processed/loading_AllWatersheds_annual.csv",
# row.names = FALSE)
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TN.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TN load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 34 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TDN.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TDN load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 105 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TP.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TP load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 24 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TDP.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TDP load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 34 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.TOC.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "TOC load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 105 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds_differences, aes(x = Month, y = load.DOC.difference)) +
geom_point(alpha = 0.7, size = 2) +
facet_wrap(vars(WS_Name), ncol = 3, scales = "free_y") +
scale_x_continuous(n.breaks = 12) +
labs(y = "DOC load difference (UTAHDWQ-WFWQC, ton/mo)", color = "Organization") +
scale_color_viridis_d(option = "magma", begin = 0.2, end = 0.7, direction = -1) +
theme(legend.position = "top")
## Warning: Removed 105 rows containing missing values (geom_point).
ggplot(loading_AllWatersheds,
aes(x = WS_Name, y = load.TN.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
facet_wrap(vars(OrganizationIdentifier)) +
labs(x = "", y = "TN loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 35 rows containing missing values (position_stack).
ggplot(loading_AllWatersheds_dmr,
aes(x = WS_Name, y = load.TKN.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "TN loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 12 rows containing missing values (position_stack).
ggplot(subset(loading_AllWatersheds, OrganizationIdentifier == "UTAHDWQ_WQX"),
aes(x = WS_Name, y = load.TDN.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "TDN loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
ggplot(loading_AllWatersheds,
aes(x = WS_Name, y = load.TP.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
facet_wrap(vars(OrganizationIdentifier)) +
labs(x = "", y = "TP loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 24 rows containing missing values (position_stack).
ggplot(loading_AllWatersheds_dmr,
aes(x = WS_Name, y = load.TP.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "TN loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
ggplot(loading_AllWatersheds,
aes(x = WS_Name, y = load.TDP.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
facet_wrap(vars(OrganizationIdentifier)) +
labs(x = "", y = "TDP loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 37 rows containing missing values (position_stack).
ggplot(loading_AllWatersheds_dmr,
aes(x = WS_Name, y = load.TDP.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "TN loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
## Warning: Removed 12 rows containing missing values (position_stack).
ggplot(subset(loading_AllWatersheds, OrganizationIdentifier == "UTAHDWQ_WQX"),
aes(x = WS_Name, y = load.TOC.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "TOC loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
ggplot(subset(loading_AllWatersheds, OrganizationIdentifier == "UTAHDWQ_WQX"),
aes(x = WS_Name, y = load.DOC.tonmo, fill = as.factor(Month))) +
geom_bar(stat = "identity") +
scale_fill_viridis_d(option = "magma") +
labs(x = "", y = "DOC loading (tons/year)", fill = "Month") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
ggplot(loading_AllWatersheds_annual_combined,
aes(x = WS_Name, y = load.TP.tonyr, fill = OrganizationIdentifier)) +
geom_bar(stat = "identity", position = position_dodge(preserve = "single")) +
scale_fill_viridis_d(option = "magma", begin = 0.2, end = 0.8) +
labs(x = "", y = "TP loading (tons/year)", fill = "Source") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
ggplot(loading_AllWatersheds_annual_combined,
aes(x = WS_Name, y = load.TN.tonyr, fill = OrganizationIdentifier)) +
geom_bar(stat = "identity", position = position_dodge(preserve = "single")) +
scale_fill_viridis_d(option = "magma", begin = 0.2, end = 0.8) +
labs(x = "", y = "TN loading (tons/year)", fill = "Source") +
theme(legend.position = "top", axis.text.x = element_text(angle = 45, hjust = 1))
TO DO: Change y axis to volumes after determining monthly lake area